Autonomous optimization of single-well and multi-well gas lift

ABSTRACT

Methods and systems are provided for controlling gas lift in at least one production well, which determine a quadratic function that relates flow rate or value of oil produced from a production well to gas flow rate for gas injected into the production well based on flow rate measurements for gas injected into the production well and corresponding flow rate measurements of oil produced from the production well. The quadratic function is used to determine an optimal flow rate for gas injected into the production well. The optimal flow rate for gas injected into the production well is used to control the flow rate of gas injected into the production well, and an oil flow rate produced from the production well is measured at the optimal gas flow rate for iterative processing to refine the quadratic function and determination of optimal gas flow rate if need be.

BACKGROUND 1. Field

The present disclosure related to artificial lift of produced fluids inone or more wells that traverse a hydrocarbon-bearing formation, and,more particularly to control of gas lift in single-well and multi-wellapplications.

2. State of the Art

An oil field can utilize an artificial lift system to lift producedfluids in a well to the surface. One such artificial lift system is agas lift system which employs a mandrel generally mounted alongproduction tubing and lowered into the well's production casing togetherwith the tubing. Gas is introduced into the annular region between thecasing and the tubing under pressure, and valves positioned along and/orwithin the mandrel allow the gas to be introduced into the fluid flowwithin the production tubing. The gas lift system helps lift producedfluids to the surface by reducing the density of the produced fluid (andthus the downhole pressure), which accelerates the movement of thefluids up the production tubing.

Oil field operators can employ a control system to manage the operationof the gas lift system in order to optimize production of fluids at thesurface. The paper by K. Rashid “Optimal Allocation Procedure forGas-Lift Optimization,” Ind. Eng. Chem. Res., 49, 2010, pp. 2286-2294,describes gas-lift optimization where lift performance is described interms of a well-head pressure and the injected gas rate. The problemthen needs to consider a network of pipelines and an offline/onlinesequence may be required.

Autonomous control and optimization of the gas lift system may becumbersome. Consideration of well-head pressures and pipelinesintroduces complexities that are difficult and costly to implement forreal-time automated applications.

SUMMARY

Methods and systems are provided for controlling gas lift in at leastone production well, which involve:

i) determining a quadratic function that relates flow rate or value ofoil produced from a production well to gas flow rate (which can beunderstood to be at specified standard conditions) for gas injected intothe production well based on flow rate measurements for gas injectedinto the production well and corresponding flow rate measurements of oilproduced from the production well;

ii) using the quadratic function to determine an optimal flow rate forgas injected into the production well; and

iii) using the optimal flow rate for gas injected into the productionwell to control the flow rate of gas injected into the production well,and measuring an oil flow rate produced from the production well at theoptimal gas flow rate for iterative processing to refine the quadraticfunction and determination of optimal gas flow rate if need be.

In embodiments, the operations of i) to iii) can be repeated, whereinfor one or more iterations of i) to iii), the quadratic function of i)is based on the measurement of oil flow rate in the previous iterationof iii).

In embodiments, the measurement of oil flow rate in the previousiteration of iii) replaces a data point used to determine the quadraticfunction in the previous iteration of iii) and preservesnon-monotonicity of the quadratic function.

In embodiments, the operations of i) to iii) can be repeated for one ormore iterations until a predetermined criterion is satisfied. Inembodiments, a shift in oil rate at the previously set optimal gasinjection rate is taken into account while preserving nonmonotonicity toupdate the quadratic function and obtain a new optimal gas injectionrate.

In embodiments, the quadratic function relates oil flow rate producedfrom the production well to injected gas flow rate for gas injected intothe production well. In other embodiments, the quadratic functionrelates oil flow rate produced from the production well to shifted gasflow rate for gas injected into the production well, wherein shifted gasflow rate is based on injected gas flow rate less a threshold gas flowrate.

In embodiments, the optimal flow rate for gas injected into theproduction well as determined in ii) can be based on coefficients of thequadratic function. For example, the optimal flow rate can be calculatedas

$- {\frac{B_{1}}{2A_{1}}.}$

In embodiments, for multi-well applications, the operations of i) can beperformed separately for a plurality of production wells to determine acorresponding plurality of quadratic functions, wherein each quadraticfunction relates oil flow rate produced from one of the plurality ofproduction wells to gas flow rate for gas injected into that oneproduction well. Operations of ii) can use the plurality of quadraticfunctions to determine optimal flow rates for gas injected into each oneof the plurality of production wells. The operations of iii) can use theplurality of optimal flow rates for gas injected into plurality of theproduction wells to control the flow rate of gas injected into theplurality of production wells, and measure oil flow rates produced fromthe plurality of production wells at the optimal flow rates of gasinjected into the plurality of production wells.

In embodiments, the optimal flow rate for gas injected into at least oneof the plurality of production wells as determined in ii) can be basedon coefficients of the quadratic function for the correspondingproduction well. For example, the optimal flow rate can be calculated as

$- {\frac{B_{1}}{2A_{1}}.}$

In embodiments, the optimal flow rate for gas injected into at least oneof the plurality of production wells as determined in ii) can be basedon a constraint. For example, the constraint can be based on a totalavailable flow rate Q_(g) of injected gas from a common gas sourceshared by the plurality of production wells. Alternatively oradditionally, the constraint can be based on a threshold rate ofinjected gas for at least one production well.

In embodiments, the optimal flow rates for gas injected into theplurality of production wells as determined in ii) can be based on asystem of nonlinear equations solved by a sequential quadraticprogramming (SQP) solver. In these embodiments, the optimal flow ratesfor gas injected into the plurality of production wells can beformulated as a nonlinear optimization problem that can be solved by asequential quadratic programming (SQP) solver.

Systems for controlling gas lift in at least one production well arealso described and claimed, which include a first flow meter formeasuring flow rate of gas injected into a production well, a controlvalve for controlling flow rate of gas injected into the productionwell, and a second flow meter for measuring flow rate of oil producedfrom the production well. A controller is operably coupled to the firstflow meter, the control valve, and the second flow meter. The controllercan be configured to perform the operations of i) to iii) to control gaslift in an optimal manner as described and claimed.

BRIEF DESCRIPTION OF THE DRAWINGS

FIGS. 1 and 2 are schematic diagrams of a single exemplary productionwell that employs a gas lift system for artificial lift.

FIG. 3 is a schematic diagram of two exemplary production wells thatemploy a gas lift system for artificial lift.

FIG. 4A is a plot of the flow rate of produced oil as a function of gasinjection flow rate for a single exemplary production well employing thegas lift system of FIG. 2. The plot shows two curves, both of which arequadratic in form. The dashed quadratic is based on field data shownwith markers and represents a standard least squares fit with equalweights for all of the data points. The solid quadratic is determinedusing the control scheme as described in the present disclosure.

FIG. 4B is a plot of the flow rate of produced oil as a function of gasinjection flow rate for three exemplary production wells employing a gaslift system similar to that shown in FIG. 3. The plot shows threecurves, all of which are quadratic in form and determined using thecontrol scheme as described in the present disclosure. The plot alsoshows field data as different markers for the three production wells.

FIGS. 5A-5E, collectively, is a flow chart that illustrates a controlscheme for gas lift operations for multiple wells based on directconstrained optimization using shifted gas flow rates.

FIG. 6 is a flow chart that illustrates adaptations to the gas liftcontrol scheme of FIGS. 5A-5E based on nonlinear optimization of gasflow rates.

FIG. 7 is a functional block diagram of an exemplary computer processingsystem.

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

The particulars shown herein are by way of example and for purposes ofillustrative discussion of the embodiments of the subject disclosureonly and are presented in the cause of providing what is believed to bethe most useful and readily understood description of the principles andconceptual aspects of the subject disclosure. In this regard, no attemptis made to show structural details in more detail than is necessary forthe fundamental understanding of the subject disclosure, the descriptiontaken with the drawings making apparent to those skilled in the art howthe several forms of the subject disclosure may be embodied in practice.Furthermore, like reference numbers and designations in the variousdrawings indicate like elements.

The systems and methods described herein operate to control operationsof a gas lift system that lifts produced fluids from one or more wells,such as those found in oil and gas fields. Such fields generally includeone or more production wells that provide access to the reservoir fluidsunderground. FIG. 1 shows an example production well 100 with a borehole102 that has been drilled into the earth. Such boreholes can be drilledto ten thousand feet or more in depth and can be steered horizontallyfor perhaps twice that distance. The production well 100 also includes acasing header 104 and casing 106, both secured into place by cement 103.Blowout preventer (BOP) 108 couples to casing header 104 and thewellhead 110, which together seal in the wellhead and enable fluids tobe extracted from the well in a safe and controlled manner.

A packer 122 and a gas lift mandrel (or tool) 114 are coupled in-linewith production tubing 112 that extends to the wellhead 110. The packer122 is configured to isolate the production zone 125 of the reservoirbelow the packer 122 from the upper part of the well. Produced fluidsflow from the production zone 125 into the flowline of the productiontubing 112. Pressurized gas is injected from the surface into theannular region between the casing 106 and the production tubing 112 andgas lift mandrel 114. The gas lift mandrel 114 has valves that areconfigured to permit injection of the gas from the annular regionbetween the casing 106 and the mandrel 114 into the flowline of themandrel 114 and the production tubing 112. In this manner, the gas isintroduced into the fluid flow within the production tubing 112 andhelps lift produced fluids to the surface by reducing the density of theproduced fluid, which accelerates the movement of the fluids up theproduction tubing 112.

FIG. 2 shows a diagram of an illustrative gas lift system incorporatedinto the production well of FIG. 1, and includes some components notshown in FIG. 1 while excluding others for clarity. Pressurized gas issupplied by a gas source (not shown) and injected into the annularregion 150 between the casing 106 and mandrel 114 under pressure via acontrol valve (or choke) and injection line flow meter 152 located atthe surface The control valve of block 152 can be operated under controlof a surface-located controller 162 to regulate the gas injectionpressure or the flow rate of the gas that is injected into the well 100during operation of the well. The injection line flow meter of block 152can be configured to measure the flow rate of the gas that is injectedinto the well 100 during operation of the well and communicate data orsignals representative of such flow rate measurements to the controller162. The communication interface(s) between the controller 162 and thevalve and injection line flow meter of block 152 can be a wiredelectrical communication link, a wireless communication link, an opticalcommunication link or other suitable communication system. The injectionline flow meter can be disposed upstream or downstream from the controlvalve. The gas lift mandrel 114 includes valves 154, 155 that allow thepressurized gas to flow from the annular region 150 into the flowline ofthe mandrel 114 and the production tubing 112, while preventing producedfluid (e.g., oil) within the flowline of the mandrel 114 and theproduction tubing 112 from flowing back out into annular region 150. Theinjected gas can mix with the produced fluid (e.g., oil) in the flowlineand reduce the density of the produced fluid. In this manner, theinjected gas can help lift produced fluids to the surface and acceleratethe movement of the produced fluids up the production tubing 112. At thesurface, the produced fluids pass through the wellhead 110 to a phaseseparator and production line flow meter 157 for output to a productionline 158. The phase separator of block 157 can be configured to separatecomponents (such as oil, water and possibly gas) of the flow of producedfluid that exits the wellhead 110. The oil component flow produced bythe phase separator flows to the production line 158. The productionline flow meter of block 157 is disposed downstream of the phaseseparator and is configured to measure the flow rate of the oil thatflows from the phase separator into and through the production line 158during operation of the well. The production line flow meter cancommunicate data or signals representative of such flow ratemeasurements to the controller 162. The communication interface betweenthe controller 162 and the production line flow meter can be a wiredelectrical communication link, a wireless communication link, an opticalcommunication link or other suitable communication system. The producedgas may be collected, compressed, and reinjected.

FIG. 3 shows a diagram of an illustrative gas lift system incorporatedinto two production wells 100 a, 100 b, which are similar to theproduction well 100 of FIGS. 1 and 2. The gas lift system includes manycomponents (including control valve and injection flow meter of block152, gas lift mandrel 114, and phase separator and production lineflowmeter of block 157) that are replicated for each one of the wells100 a, 100 b and operate under control of a common controller 162. Notethat the produced oil that flows through the production lines 158 of thetwo wells 100 a, 100 b can be combined for output to a pipeline asshown. Furthermore, the pressurized gas that is supplied to the controlvalves of block 152 of the two wells 100 a, 100 b can be provided by acommon gas source (not shown). In this case, because the pressurized gasused for the gas lift is shared between multiple wells, the maximum gasflow into the wells can be constrained by the maximum output of thecommon gas source. Such constraints can be limiting when a number ofwells share the output of a common gas source.

Gas-Lift Principles

When a production well is unable to lift produced fluid (liquidhydrocarbon) over the entire vertical height of the well (H), thereservoir pressure is lower than the wellhead pressure (P₀) plus thepressure due to the column height. Artificial lift is required if

$\begin{matrix}{{{P_{0} + {\int\limits_{- H}^{0}{\rho_{o}{gdz}}}} > P_{f}},} & {{Eqn}.\mspace{14mu} (1)}\end{matrix}$

where P_(f) is the formation pressure, z is the vertical heightreferenced to the wellhead, ρ_(o) is the density of the produced fluid(liquid hydrocarbon), and g is the acceleration due to gravity. Notethat the integration is with respect to a vertical coordinate z.

Alternatively, an equilibrium height (H_(e)) may be obtained by solving

$\begin{matrix}{{{P_{0} + {\int\limits_{- H_{e}}^{0}{\rho_{o}{gdz}}}} = P_{f}},} & {{Eqn}.\mspace{14mu} (2)}\end{matrix}$

for H_(e). Obviously given that ρ_(o)>0, H_(e)<H when artificial lift isnecessary. In general, ρ_(o) varies with pressure p and temperature T,but in comparison to gas density, ρ_(g), this may be a minor variation.

The principle in gas lift is to inject gas at a particular flow rate(q_(g)), preferably optimally, in order to achieve as high an oilproduction flow rate (q_(o)) as possible by reducing the fluid mixturedensity of the produced fluid and gas lower than q_(o) so that H_(e)>H.In embodiments, the flow rate q_(g) for a particular well can bedetermined from flow rate measurement(s) of a surface-located gasinjection line flowmeter (e.g., block 152 of FIGS. 2 and 3), reported atstandard conditions. Such standard conditions as per IUPAC, refer to atemperature and pressure of 0° C. and 100 kPa respectively. The oilproduction flow rate (q_(o)) can be determined from the flow ratemeasurement(s) of one or more a surface-located production lineflowmeters (e.g., block 157 of FIGS. 2 and 3). Note that a low gasinjection flow rate does not reduce the produced fluid densitysufficiently to achieve the optimal rate, and too high a gas injectionrate reduces the liquid fraction in the produced fluid so much that theproduced oil flow rate suffers. More importantly, since H_(e)<H, q_(g)has to be above a threshold flow rate (θ_(g)) for the oil productionflow rate q_(o)>0. At a given reservoir state, there then exists anoptimal flow rate of injected gas (denoted q_(g,O)) which can bespecified as

$\begin{matrix}{{\frac{d{q_{o}\left( q_{g,O} \right)}}{dq_{g}} = 0},} & {{Eqn}.\mspace{11mu} (3)}\end{matrix}$

at the current reservoir condition. This optimal flow rate q_(g,O) isgreater than θ_(g), and is always positive. A constraint on the totalavailable flow rate of injected gas can be given and denoted Q_(g). Fora single well, this constraint is not an issue. However, in multi-wellapplications, this constraint can be of issue and accounted for asdescribed below.

Usually, the reservoir time scales are much larger than the wellboreequilibration time, and for the purpose of this disclosure, thereservoir can be deemed to be pseudo-static. This separation of timescales allows the control schemes as described herein to construct aformal solution to the gas-lift optimization problem. In thesesolutions, no well is affected by an adjusted production in other wellsover a time-scale of stabilization within the well.

Advantageously, the control schemes disclosed herein do not considerwell-head pressures. Specifically, by configuring a surface-locatedcontrol valve to set the flow rate of injected gas at the wellhead, theflow rate of the oil component of the produced fluid at the surface canbe considered to be a function of the flow rate of the injected gas. Inembodiments, the flow rate of the oil component of the produced fluid atthe surface can be considered to be a function of the flow rate of theinjected gas at standard conditions. Thus, considerations of pipelinesare not needed. Furthermore, the control schemes do not require detailedmodeling of either the wellbore or the reservoir, and relies onmeasurements of flow rates at the surface and some basic principles ofthe gas-lift mechanism. Since the control schemes do not rely upondetailed simulations, the control schemes can be easily implemented forreal-time automated applications. Furthermore, the control schemes candetermine optimal set points for controlling gas injection in multiplewells independently from one another with due inclusion of limitationsimposed by constraints.

Single Well Optimization

Before constructing a general algorithmic solution, it is useful toillustrate an observation-based autonomous control scheme through asingle-well example. For many reasons that will become obvious below, itis convenient to represent the stabilized oil flow rate at the surfacewith respect to the flow rate of injected gas as a quadratic functionwhere the coefficients of the quadratic function vary with time albeitslowly in comparison to the optimization and stabilization time. Therelation of the stabilized oil flow rate with respect to the flow rateof injected gas may be represented by a higher order polynomialfunction, but a quadratic function is the order necessary and sufficientto obtain a unique optimal flow rate from Eqn. 3. A quadratic functionalso affords simple algebraic expressions for fast computations forautonomous control systems. Furthermore, the quadratic function can belocalized or solved from most recent measurements, which allows thequadratic function to be pseudo-static and it is also expected to beaccurate near the operating point.

It is to be understood that the principles discussed here may begeneralized to higher order functions of expressing q_(o) with respectto q_(g), by adding more stabilized rates to fit a polynomial, or amin-max polynomial, or a best-fit polynomial function. But a higherorder polynomial function is undesirable for many reasons: (i) a numberof (q_(g), q_(o)) data pairs need to be obtained; (ii) with increasingnumber of data pairs, the product of number of pairs to wellborestabilization time may approach a non-negligible fraction ofinterference time between wells; and (iii) expressions become needlesslycomplex with the possibility of multiple extrema.

Note that a quadratic function has the benefit that only three datapoints are needed at any stage of a decision step. Thus, for one well,with integer i∈1,3 representing the three data points, the data pointsfor the optimal flow rate q_(g,O) can be given as

g _(o1i) =A ₁ q _(g1i) ² +B ₁ q _(g1i) +C ₁,  Eqn. (4)

where it is assumed that the A₁, B₁, and C₁ are updated when a new pointis added to the flow rates and an old one dropped. The subscript 1 forA₁, B₁, and C₁ is an index that refers to the well, and these values arederived from the most recent data. With regard to the notation used inEqn. (4), the first subscript o for the flow rate q_(o) denotes the oilphase, the second subscript 1 for the flow rate q_(o) represents thewell, and the third subscript i for the flow rate q_(o) denotes the datapoint. Similarly, the first subscript g for the flow rate q_(g) denotesthe gas phase, the second subscript 1 for the flow rate q_(g) representsthe well, and the third subscript i for the flow rate q_(g) denotes thedata point. For instances where the third subscript is absent impliesonly the well, and instances where both the second and third subscriptsare absent refers to the phase in a general context.

Given three data points for the flow rates, the coefficients are

$\begin{matrix}{\mspace{79mu} {{A_{1} = {\frac{1}{q_{g12} - q_{g13}}\left\lbrack {\frac{q_{o\; 11} - q_{o\; 12}}{q_{g11} - q_{g12}} - \frac{q_{o\; 11} - q_{o\; 13}}{q_{g11} - q_{g13}}} \right\rbrack}},}} & {{Eqn}.\mspace{11mu} (5)} \\{{B_{1} = {\frac{q_{o\; 11} - q_{012}}{q_{g11} - q_{g12}} - \frac{\begin{matrix}{q_{g11} + q_{g12}} \\\begin{bmatrix}{{q_{g11}\left( {q_{o\; 13} - q_{012}} \right)} +} \\{{q_{g13}\left( {q_{o\; 12} - q_{011}} \right)} + {q_{g12}\left( {q_{o\; 11} - q_{013}} \right)}}\end{bmatrix}\end{matrix}}{\left( {q_{g12} - q_{g13}} \right)\left( {q_{g11} - q_{g12}} \right)\left( {q_{g11} - q_{g13}} \right)}}},\mspace{20mu} {and}} & {{Eqn}.\mspace{11mu} (6)} \\{\mspace{79mu} {C_{1} = {q_{o\; 11} - {A_{1}q_{g11}^{2}} - {B_{1}{q_{g11}.}}}}} & {{Eqn}.\mspace{11mu} (7)}\end{matrix}$

Knowing A₁ and B₁ from the three settings of q_(g) is now enough to setthe optimal rate of q_(g) (for A₁≠0) as

$\begin{matrix}{{q_{{g\; 1},O} = {- \frac{B_{1}}{2A_{1}}}}.} & {{Eqn}.\mspace{11mu} (8)}\end{matrix}$

Eqn. (5) shows that A₁<0, thereby assuring upward convexity, if thethree points are non-monotonic for q_(o) with increasing q_(g). If A₁and B₁ are both zero, or close to it, for the ratio to be meaningless,the flow rate q_(o1) is a constant, and q_(g1) does not affect q_(o1)suggesting that all of the three rates are acceptable and are close tooptimal. To reemphasize, in the absence of measurement errors, q_(o) isa convex-upward function, and therefore A₁<0.

If q_(g1,O) is substantially different from the first three rates, anoption to limit the next setting may be imposed so that q_(g1,O) doesnot change from min(q_(g1)) or max(q_(g1)) by a preset threshold value.Now, upon setting q_(g1,O), or a change-limited version of it, a newdata point is obtained. Data point i=1 (the earliest) may be dropped,and the others moved down in index with q_(g1,O) becoming q_(g13).Another option is to search the three points and discard the point thatstill preserves non-monotonicity in q_(o) vs. q_(g). This approach keepsthe new q_(g1,O) within bounds. The process may then be continued, untilan acceptable convergence to q_(g1,O) is found. Once thus determined,the flow rate q_(g1) needs to be periodically set to a different valueaway from the current q_(g1,O) (possibly one of q_(g11) and q_(g12)) forupdating the parameters for A₁, B₁, and C₁ of the quadratic equation andthe solution for the new q_(g1,O). The timescale for such periodicupdate to the new flow rate q_(g1) can be selected to be greater thanthe characteristic interwell interference time for the field in order toallow the system to equilibrate. Note that if q_(g1,O)>Q_(g), the flowrate q_(g1,O) can be limited to Q_(g). In this case, the oil productionrate is limited by gas supply. Note that the optimal gas injectionrate(s) as determined from the methodology described herein need not bethe most optimal injection rate(s) but can provide a solution that isclose to the most optimal injection rate within an acceptable tolerancerange.

The description above demonstrates how single-well gas lift optimizationmay be carried out purely by observation. No reservoir knowledge otherthan to know the separation of well-bore time scale from that of theinterwell interference time is necessary.

Two-Well System

A similar control scheme can be applied for simultaneous gas-liftoptimization of two wells. Because the time scales of interference aremuch larger than the time scale of well response to gas injection, thetwo-dimensional optimization problem can be disentangled into twoone-dimensional optimizations subject to imposed constraints. Commonpipelines for the handling of oil can be assumed to be large enough thatthere is no limit on the produced oil flow rate. Because the controlscheme for gas rate optimization is based on flow control, it isreasonable to assume that any common pressure link to gas lines isinconsequential. For example, an extra drawdown in one gas line willreduce the gas supply pressure, but the second line's control valveshould open correspondingly to maintain the requisite gas injection flowrate.

It is evident that as long as there are no bounds imposed on the totalflow rate of oil Q_(o) from the two wells and no bounds imposed on thetotal flow rate of injected gas into the two wells given as Q_(g), thecontrol of the gas injection flow rates q_(g1,O) and q_(g2,O) shouldyield the maximum q_(o) as long as q_(g1,O)≥θ_(g1) and q_(g1,O)≥θ_(g2),where θ_(g1) or θ_(g2) are threshold gas injection flow rates for thetwo wells. If not, the injection gas flow rate q_(g1) or q_(g2) can beset to the corresponding threshold flow rate θ_(g1) or θ_(g2) wheneverthe inequality is not satisfied, thereby accepting a lower q_(o) thanindicated without constraints (for example, a negative q_(g1) or q_(g2)is not possible).

The control scheme gets more complicated if a limit is set on Q_(g). Inthis case, it is assumed that

Q _(g)>θ_(g1)+θ_(g2).  Eqn. (9)

If Q_(g) is greater than [max(q_(g1,O), θ_(g1))+max(q_(g2,O), θ_(g2))],the optimal points are not constrained, and the wells can operateindependently just as in a single well. Suppose Q_(g) is less than[max(q_(g1,O), θ_(g1))+max(q_(g2,O), θ_(g2))] as determined by theindependent well assumption. Ignoring the threshold requirement for themoment, with the bound Q_(g) specified, the optimum operation is onesuch that q_(g1)+q_(g2)=Q_(g), but Q_(o)=q_(o1)+q_(o2) is maximized. Insuch cases, the control scheme can first determine A₁, A₂, B₁, B₂, C₁and C₂, either by operating one well at a time, or choosing injectiongas flow rates for the two wells such that their sum ≤Q_(g).

Optimizing q_(o) with respect to one of the injection gas flow ratessubject to q_(g1)+q_(g2)=Q_(g) gives the optimal rate

$\begin{matrix}{{q_{{g\; 1},O} = {- \frac{{- \left( {B_{1} - B_{2}} \right)} + {2A_{2}Q_{g}}}{2\left( {A_{1} + A_{2}} \right)}}},} & {{Eqn}.\mspace{11mu} (10)}\end{matrix}$

and, by symmetry,

$\begin{matrix}{{q_{{g\; 2},O} = {- \frac{{- \left( {B_{2} - B_{1}} \right)} + {2A_{1}Q_{g}}}{2\left( {A_{1} + A_{2}} \right)}}}.} & {{Eqn}.\mspace{11mu} (11)}\end{matrix}$

When thresholds are not violated, the gas injection flow rates for thetwo wells can be set to these values q_(g1,O), q_(g2,O). And theresulting oil production flow rates can be measured and the parametersA₁, A₂, B₁, B₂ for the two quadratic functions may be updated. Insubsequent steps, q_(g1) and q_(g2) can be calculated as though there isno bound, and subsequently check whether the sum q_(g1)+q_(g2) exceedsQ_(g). If so, the optimal gas injection flow rates q_(g1,O) and q_(g2,O)can be calculated according to Eqns. (10) and (11), respectively, butwith the updated quadratic functions for q_(g)→{circumflex over(q)}_(o)(q_(g)), where {circumflex over ( )} denotes the quadraticfunction.

When the thresholds are violated by Eqns. (10) or (11), the flow rate ofinjected gas for that well can be set to the corresponding thresholdrate (θ_(g1) or θ_(g2)). Without loss of generality, sayq_(g2,O)<θ_(g2), the operating condition for the second well can be setto

$\begin{matrix}{{q_{{g2},O} = {\max \left( {{- \frac{{- \left( {B_{2} - B_{1}} \right)} + {2A_{1}Q_{g}}}{2\left( {A_{1} + A_{2}} \right)}},\theta_{g2}} \right)}},} & {{Eqn}.\mspace{11mu} (12)} \\{q_{{g1},O} = {{\min \left( {{Q_{q} - q_{{g2},O}},{- \frac{B_{1}}{2A_{1}}}} \right)}.}} & {{Eqn}.\mspace{11mu} (13)}\end{matrix}$

Thus, when both the threshold and limit on total injection gas flow rateare imposed, given {circumflex over (q)}_(o)(q_(g)) for the two wells,Eqns. (10) and (11) can be used to compute and set the two gas injectionflow rates if q_(g1,O)+q_(g2,O) is less than Q_(g). Otherwise, whenq_(g1,O)+q_(g2,O) exceeds Q_(g), and if one of them from Eqns. (10) and(11) violates the bound of the threshold θ_(gi) (both cannot violatesince the sum of the thresholds <Q_(g)), then the gas injection flowrate for that well can be set to the corresponding threshold θ_(g1) orθ_(g2), and the gas injection flow rate for the other well can becalculated and set according to Eqn. (13) where the index for the onewell set to the threshold rate is indicated as well 2 in Eqn. (13). Uponstabilization, the quadratic functions for the two wells can be updatedand the process can be repeated.

Multiple Wells

For applications involving gas lift of multiple wells, a notation isemployed where a single subscript index i for the flow rates andquadratic coefficients A, B and C refers to well i. As described in thetwo-well problem, from three different settings of flow rates q_(g) fora given well i, the flow rate of produced oil can be given as

q _(oi) =A _(i) q _(gi) ² +B _(i) q _(gi) +C _(i).  Eqn. (14)

Note that A_(i)<0. Optimizing the wells, without any restraint on theamount of available gas, gives

$\begin{matrix}{q_{{gi},O} = {{\max \left( {{- \frac{B_{i}}{2A}},\theta_{gi}} \right)}.}} & {{Eqn}.\mspace{11mu} (15)}\end{matrix}$

Obviously, the thresholds θ_(gi) can be set so that at θ_(gi), q_(gi)≥0(usually q_(0i)=0). After the thresholds are met, it is possible forΣ_(i)q_(gi,O)>Q_(g). If not, all of the gas injection flow rates can beset to q_(g1,O) and the process continues. If the constraint in Q_(g) isviolated, q_(g1) and q_(g2) can be set such that the sum of q_(g1) andq_(g2) is equal to Q_(g). Now, with the equality constraint, writing tominimize subtraction errors, it is preferable (and by no meansnecessary) to choose the gas injection flow rate q_(gN) for the lastwell N to be the largest expected one as follows

$\begin{matrix}{q_{qN} = {Q_{g} - {\sum\limits_{= 1}^{N}{q_{gi}.}}}} & {{Eqn}.\mspace{11mu} (16)}\end{matrix}$

This allows for calculation of the total flow rate of oil production asfollows

$\begin{matrix}{Q_{o} = {\left( {\sum\limits_{j = 1}^{N - 1}\left( {{A_{j}q_{gj}^{2}} + {B_{j}q_{gj}} + C_{j}}\  \right)} \right) + {q_{oN}.}}} & {{Eqn}.\mspace{11mu} (17)}\end{matrix}$

For the moment, the thresholds can be ignored but with a requirementthat the sum of gas injection flow rates does not exceed Q_(g). If thewells are numbered from 1 to N where N>2, the gas injection flow rate ofwell Nin Eqn. (17) can be replaced as follows

$\begin{matrix}{Q_{o} = {\left( {\sum\limits_{j = 1}^{N - 1}\left( {{A_{j}q_{gj}^{2}} + {B_{j}q_{gj}} + C_{j}}\  \right)} \right) + {\left( {{A_{N}q_{gN}^{2}} + {B_{N}q_{gN}} + C_{N}} \right).}}} & {{Eqn}.\mspace{11mu} (18)}\end{matrix}$

Maximizing Q_(o) with respect to q_(gi) yields N−1 relationships asfollows

$\begin{matrix}{{{{2A_{i}q_{{gi},O}} + B_{i} - {2{A_{N}\left( {Q_{g} - {\sum\limits_{j = 1}^{N - 1}q_{{gi},O}}} \right)}} - B_{N}} = 0}{{{{for}\mspace{20mu} i} = 1},2,\ldots \;,{N - 1.}}} & {{Eqn}.\mspace{11mu} (19)}\end{matrix}$

The linear system of equations of (19) provides solutions for q_(gi,O)for i arranging from 1 to N−1. And from Eqn. (16), the optimal q_(gN)for the last well N is also known. In matrix form, Eqn. (19) may berepresented as

Vu=w,  Eqn. (20)

where the elements of V are

v _(ij)=2A _(i)δ_(ij)+2A _(N),  Eqn. (21)

where δ_(ij) is the Kroneckar delta, andwhere the elements of u are

u _(i) =q _(gi,O), and  Eqn. (22a)

where the elements of w are

w _(i) =B _(N) −B _(i)+2A _(N) Q _(g).  Eqn. (22b)

The q_(gi,O) of Eqn. (22a) can be obtained by solving the linear systemof equations of (20). Naturally, the sum of all the flow rates q_(gi,O)would be Q_(g). A check is however necessary for any violation of thethreshold constraint imposed by θ_(gi). For example, if M of the flowrates q_(gi,O) violate this constraint, it is quite possible that thepreviously indexed Nth flow rate could be one of the M. If it does notsatisfy the lower bound θ_(gN) constraint, the indices for the violatingwells run from (N−M+1) to N. In this case, the flow rates q_(gi) for theviolating wells (N−M+1) to N can be set to θ_(gi). If the N^(th) wellsatisfies the constraint, then the index for this well is shifted to(N−M) and again set the flow rates for the violating wells (N−M+1) to Nto θ_(gi). The previous Q_(g) can be used to define the flow rate Q_(g0)and Q_(g) can be updated to be

$\begin{matrix}{Q_{g} = {Q_{g0} - {\sum\limits_{i = {N - M + 1}}^{N}{\theta_{gi}.}}}} & {{Eqn}.\mspace{11mu} (23)}\end{matrix}$

Now the optimization is reduced to N−M wells, but satisfying theconstraint that the total gas injection rate must be at most Q_(g0).Since all of the wells greater than an index of N−M have had their flowrates increased from optimum setting in order to meet the threshold. theold optimal positions for wells indexed from 1 to N−M will result in atotal rate exceeding Q_(g0). Therefore, it is necessary that theequality constraint for wells indexed from 1 to N−M be satisfied. Thus,we may again eliminate well with index N−M because

$\begin{matrix}{{Q_{g} = {\sum\limits_{i = 1}^{N - M}q_{gi}}},} & {{Eqn}.\mspace{11mu} (24)}\end{matrix}$

where Q_(g) satisfies Eqn. (23).Eqn. (19) can be rewritten as

$\begin{matrix}{{{{2A_{i}q_{{gi},O}} + B_{i} - {2{A_{N - M}\left( {Q_{g} - {\sum\limits_{j = 1}^{N - M - 1}q_{{gi},O}}} \right)}} - B_{N - M}} = 0}{{{{for}\mspace{14mu} i} = 1},2,\ldots \;,{N - M - 1.}}} & {{Eqn}.\mspace{11mu} (25)}\end{matrix}$

Note that Eqns. (21), (22(a)) and (22(b) are still applicable along withEqn. (20) for solving for the N−M−1 unknowns for the flow rates q_(gi)with the (N−M)^(th) q_(gi) obtained from Eqn. (24) and A_(N-M) andB_(N-M) replacing A_(N) and B_(N). A check is again made to see whetherany of the θ_(gi) is violated by the N−M values for q_(gi). For any ofthe violated q_(gi), the gas injection flow rate can be set to thecorresponding threshold θ_(gi) and the process can be repeated until allthe optimum values meet the constraint as well.

With all of the flow rates set for q_(gi), new values of q_(oi) can bemeasured and the quadratic functions for the wells updated. As part ofupdating the quadratic functions, one of the three datapoints thatdefine the quadratic function for a given well can be dropped andreplaced by the new data point given by the q_(gi) value and themeasured q_(oi) value. Preferably, the datapoint that is dropped is theoldest datapoint that still preserves non-monotonicity of the q_(o) vs.q_(g) relation. New values for the A_(i), B_(i) and C_(i) coefficientsof the quadratic functions for the wells can also be calculated. Theprocess may be continued ad infinitum so that all wells operate asoptimally as possible over time.

Note that the multi-well control scheme does not require reservoircharacterization and model building. The process may be automated.

Single Well Example

The control scheme can be illustrated with an example using a singlewell. This example assumes that there is no constraint on Q_(g) for thesingle well. Based on representative field data, an example oilproduction flow rates vs. gas injection flow rates is shown in in FIG.4A. The plot shows two curves, both of which are quadratic. The dashedquadratic uses all of the data shown with markers and represents astandard least squares fit with equal weights for all of the points. Thesolid curve uses three data points represented by filled markers, and inthis case are at approximately q_(g1)={99.11, 169.9, 722.1} ML d⁻¹. Thequadratic curve with oil rate (in m³ d⁻¹) is of the form

q _(o1)=−0.0006286q _(g1) ²+0.6183q _(g1)+135.7.  Eqn. (26)

Eqn. (8) would set

q _(g1,O)=491.8 ML d⁻¹.  Eqn. (27)

Assuming that the separation of time scales discussed earlier stillholds, the new oil rate that would be obtained is 268.2 m³ d⁻¹. Droppingthe first point at a q_(g) of 99.11 ML d⁻¹, the new (rounded) gas ratesused for the quadratic would be q_(g)i={491.8, 169.9, 722.1}MLd⁻¹ with

q _(o1)=−0.0003653q _(g1) ²+0.2835_(g1)+168.0.  Eqn. (28)

The revised optimal rate with the new quadratic is

q _(g1,O)=524.9 ML d⁻¹,  Eqn. (29)

for which the oil production rate as per the data of FIG. 7 would be266.7 m³ d⁻¹. Now. in order to preserve non-monotonicity in the q_(o)vs. q_(g) relation, the gas rate data point (722 MLd⁻¹) can be droppedand replaced with the new gas rate and corresponding oil rate. The gasrates then are {169.9, 524.9, 491.8} ML d⁻¹. The revised quadratic is

q _(o1)−0.0005286q _(g1) ²+0.4915q _(g1)+154.3.  Eqn. (30)

The corresponding setting for optimal gas rate is

q _(g1,O)=464.9 ML d⁻¹.  Eqn. (31)

The oil rate per FIG. 7 is 269.0 m³d⁻¹, which is demonstrably thehighest oil rate achieved so far. The new gas rates become {169.9,521.8, 464.9} MLd⁻¹. A revised quadratic then is

q _(o1)=−0.0005522q _(g1) ²+0.5079q _(g1)+152.2,  Eqn. (32)

which leads us to

q _(g1,O)=459.9 ML d⁻¹,  Eqn. (33)

for which the oil production (see FIG. 7 and interpolate) is 269.2 m³d⁻¹nearly the same as the previous step. The gas rates in the sequence are{169.9 521.8, 459.9} MLd⁻¹ with the quadratic becoming

q _(o1)=−0.0005620q _(g1) ²+0.5147q _(g1)+151.4,  Eqn. (34)

so that an updated gas rate setting becomes

q _(g1,O)=457.9 ML d⁻¹,  Eqn. (35)

with the corresponding oil rate being 269.2 m³ d⁻¹, same as the previousvalue within truncation errors. Thus, a satisfactory result is obtainedby a totally autonomous operation. The true maximum is about 270.2m³d⁻¹. The reason for the difference between values is that one of thegas rates is spaced far from the maximum.

Given the minor variation in oil production rate, as an alternative tothe above iterative sequence, at some point one may choose to stopfurther optimization, other than a periodic shift in gas injection ratein order to test the migration of the optimal operating point.Alternatively, since we appear to be always on the right side’ (higherq_(g)) of the maximum in the updates, the point corresponding to the gasrate of 169.9 ML d⁻¹ is not replaced. Thus, we have an anchoring problemthat prevents accelerated convergence. So periodically. it is advisableto migrate a distance in q_(g) (for example, say a maximum of 10% of theprevious operating point in a direction towards an anchored point.

As an example. let us choose the next gas rate to be 10% lower than theprevious setting (i.e., towards the anchored point). The rate is then413.9 ML d⁻¹ for which the oil rate is 269.6 m³ d⁻¹. To preservenon-monotonicity. the anchor point stays put, but with a new quadraticof

q _(o1)=−0.0006169q _(g1) ²+0.5529q _(g1)+146.5,  Eqn. (36)

setting the new optimal gas injection rate of

q _(g1,O)=448.1 ML d⁻¹,  Eqn. (37)

for an oil rate of 269.6 m³ d⁻¹. Since this suggests that the optimalrate is between 413.9 ML d⁻¹ and 448.1 ML d⁻¹, the prior value ofq_(g1)=169.9 ML d⁻¹ may now be discarded. Convergence to the truemaximum is reached rapidly. The new optimal gas rate with the updatedquadratic is 431.1 ML d⁻¹. with an oil production rate of 270.1 m³ d⁻¹.This is very close to the true maximum of 270.3 m³ d⁻¹ and is achievedin just three steps.

Multi-Well Example

In order to illustrate multi-well applications, consider the case ofthree wells with each well one producing oil only at a non-zero positivegas injection flow rate. The ‘data’ points and the first three ratesused for the quadratic are shown in FIG. 4B. For the purpose ofillustration, the constraints θ_(gi) are set to 30 ML d⁻¹ and Q_(g) isset to 300 ML d⁻¹. This prevents the gas flow rates being set so thatmaximum oil rate is not possible in all three wells. The minimum rateensures non-negative oil production in all wells. For the startingquadratic functions,

q _(g1)={48.14,84.95,198.2} ML d⁻¹,  Eqn. (38)

q _(g2)={48.14 84.95,198.2} ML d⁻¹, and  Eqn. (39)

q _(g3)={45.31,90.61,283.2} ML d⁻¹,  Eqn. (40)

so that

q _(o1)=−0.0001439q _(g1) ²+0.4075q _(g1)+13.93,  Eqn. (41)

q _(o2)=−0.0001748q _(g2) ²+0.5047q _(g2)−12.3, and  Eqn. (42)

g _(o3)=−0.0003389q _(g3) ²+0.1250q _(g3)−0.9936.  Eqn. (43)

Here, the combined optimal points exceed Q_(g). Using the multi-welloptimization with the quadratics and the constraint Q_(g), the optimalrates are calculated as

q _(g1,O)=113.5 ML d⁻¹,  Eqn. (44)

q _(g2,O)=121.3 ML d⁻¹, and  Eqn. (45)

q _(g3,O)=65.2 ML d⁻¹.  Eqn. (46)

Through linear interpolation between the data points shown in FIG. 8,the total oil production rate Q_(o) can be calculated as 69.97 m³ d⁻¹.Constructing updated quadratics with these new rates in whichnon-monotonicity is retained, the newly optimized gas rates are given as

q _(g1,O)=116.7 ML d⁻¹,  Eqn. (47)

q _(g2,O)=114.0 ML d⁻¹, and  Eqn. (48)

q _(g3,O)=69.3 ML d⁻¹.  Eqn. (49)

The total oil production rate Q_(o) for the newly chosen operatingpoints can be calculated as 70.24 m³ d⁻¹. The difference in the totaloil production rate Q_(o) though improved is within the margin of errorof interpolation errors, and may be considered to be close to optimaljust after the first step.

The operating point may be adjusted further when a meaningful change isdetected in lift performance. Such meaningful change can be inferredfrom a noticeable change in oil production rate (for example, 5% ormore). A new quadratic may be constructed with the updated gas flowrate, and this will automatically move the operation to a new optimalpoint. Thus, a continuous tracking of the optimal operating point isachieved.

Initialization

The automated control scheme employs three data points of the oilproduction flow rate q_(oi) as related to the gas injection flow rateq_(gi) for each well in order to obtain the initial quadratic functions.Since this is only an initialization process, and is refined with data,it is reasonable to start with q_(oi)=0 for q_(gi)=0 for each well i.The process can then start with an initial rate q_(g) and then increaseit by a factor f>1. If upon increasing it to f q_(g), and q_(o) is notreduced, an enhanced f factor can be applied, dropping the earliest rate(here it is 0), until three rates are available that result in anon-monotonic quadratic function {circumflex over (q)}_(o)(q_(g)).

Why Quadratic?

The quadratic function is the simplest form of a polynomial that allowsfor a maximum in the relation of the oil production flow rate q_(o) tothe gas injection flow rate q_(g). Furthermore, an expression for theposition of the maximum relies on the most current data of three points,and thus avoids the problem of the slow drift due to reservoir pressurechange expected over time. The separation of well performance from thatof the reservoir makes it unnecessary to utilizereservoir-model-based-optimization that is usually fraught withuncertainty.

Direct Constrained Optimization

The previous sections illustrate how, with a Q_(g) constraint, a forwardcalculation based on Eqn. (19) provides a system of linear equations forq_(gi). With θ_(gi) imposed, the problem can be reduced to sequentialreduction of constraints as done through Eqn. (19). An alternativeapproach is based on shifted gas flow rates by recognizing the following

$\begin{matrix}{{{\sum\limits_{i = 1}^{N}\theta_{gi}} < Q_{g}},{and}} & {{Eqn}.\mspace{11mu} (50)} \\{{{\overset{\hat{}}{q}}_{oi}\left( \theta_{gi} \right)} \geq 0.} & {{Eqn}.\mspace{11mu} (51)}\end{matrix}$

A shifted gas flow rate can be defined as

s _(gi) =q _(gi)−θ_(gi).  Eqn. (52)

The quadratic functions {circumflex over (q)}_(oi) can be redefined interms of s_(gi) so that

{tilde over (q)} _(oi)(s _(gi))={circumflex over (q)} _(oi)(q_(gi)).  Eqn. (53)

We also note that the true optimal operating point (denoted by anasterisk) should satisfy

q _(gi,O)*>θ_(gi) or s _(gi,O)*>0.  Eqn. (54)

Now, the function {tilde over (q)}_(oi)(s_(gi)) is strictly convexupward. With ′ indicating differentiation “with respect to” argument,{tilde over (q)}_(oi)′(s_(gi))<0. Also {tilde over (q)}_(oi)′(0)>0.

At any stage of iteration, three points s_(gi1), s_(gi2) and s_(gi3) canbe collected through which is formed a quadratic function. The threepoints s_(gi1), s_(gi2) and s_(gi3) are chosen so that oil rates withrespect to s_(gi) are non-monotonic. This implies that

$\begin{matrix}{{\frac{{q_{oi}\left( s_{{gi}\; 2} \right)} - {q_{oi}\left( s_{{gi}\; 1} \right)}}{s_{{gi}\; 2} - s_{{gi}\; 1}} > 0},{and}} & {{Eqn}.\mspace{14mu} (55)} \\{\frac{{q_{oi}\left( s_{{gi}\; 3} \right)} - {q_{oi}\left( s_{{gi}\; 3} \right)}}{s_{{gi}\; 3} - s_{{gi}\; 2}} < 0.} & {{Eqn}.\mspace{14mu} (56)}\end{matrix}$

Let the quadratic passing through the three points s_(gi1), s_(gi2) ands_(gi3) be

{tilde over (q)} _(oi) =Ã _(i) s _(gi) ² +{tilde over (B)} _(i) s _(gi)+{tilde over (C)} _(i).  Eqn. (57)

We have already shown that Ã_(i)<0. Therefore, the quadratic is alsoconvex upward. Given the mean value theorem, and the convexity of q_(oi)with respect to s_(gi), it is necessary that

$\begin{matrix}{{{{\overset{\sim}{q}}_{oi}^{\prime}\left( s_{{gi}\; 1} \right)} > \frac{{q_{oi}\left( s_{{gi}\; 2} \right)} - {q_{oi}\left( s_{{gi}\; 1} \right)}}{s_{{gi}\; 2} - s_{{gi}\; 1}} > 0},{and}} & {{Eqn}.\mspace{14mu} (58)} \\{{{\overset{\sim}{q}}_{oi}^{\prime}\left( s_{{gi}\; 3} \right)} > \frac{{q_{oi}\left( s_{{gi}\; 3} \right)} - {q_{oi}\left( s_{{gi}\; 2} \right)}}{s_{{gi}\; 3} - s_{{gi}\; 2}} < 0.} & {{Eqn}.\mspace{14mu} (59)}\end{matrix}$

Thus, the new operating point both with respect to the quadratic

$\left( {- \frac{{\overset{\sim}{B}}_{i}}{2{\overset{\sim}{A}}_{i}}} \right)$

and the true optimum determined by {tilde over (q)}_(oi)′(s_(gi))=0 issuch that

s _(gi1) <s _(gi,O) <s _(gi3).  Eqn. (60)

Since every update in iteration ensures non-monotonicity, the newlychosen point is bounded by points 1 and 3. Thus, when Q_(g) is notconstraining the optimal gas rates, the solutions are well-founded andlie between s_(gi1) and s_(gi3) and satisfy the θ_(gi) bounds.

Let us now construct optimal operating points for multiple wells with aQ_(g) constraint. Given Q_(g) and θ_(gi), noting Eqn. (50), a constraint{tilde over (Q)}_(g) can be defined as

$\begin{matrix}{{\overset{\sim}{Q}}_{g} = {{Q_{g} - {\sum\limits_{i = 1}^{N}\theta_{gi}}} > 0.}} & {{Eqn}.\mspace{14mu} (61)}\end{matrix}$

The shifted gas flow rate for the last well N can be defined as

$\begin{matrix}{s_{gN} = {{\overset{\sim}{Q}}_{g} - {\sum\limits_{i - 1}^{N - 1}{s_{gi}.}}}} & {{Eqn}.\mspace{14mu} (62)}\end{matrix}$

To determine the quadratics, the process starts with all s_(gi1)>0. Withthe quadratics obtained with s_(gi) as the independent variable, thetotal production rate (see Eqn. (19)) when maximized gives

$\begin{matrix}{{{{2{\overset{\sim}{A}}_{i}s_{{gi},O}} + {\overset{\sim}{B}}_{i} - {2{A_{N}\left( {{\overset{\sim}{Q}}_{g} - {\sum\limits_{j = 1}^{N - 1}s_{{gj},O}}} \right)}} - {\overset{\sim}{B}}_{N}} = {0\mspace{14mu} {for}}}\text{}{{i = 1},2,\ldots \mspace{14mu},{N - 1},}} & {{Eqn}.\mspace{14mu} (63)}\end{matrix}$

in which all Ã_(i)<0 and all {tilde over (B)}_(i)>0. In order to satisfyθ_(gi) bound, s_(gi,O) must be ≥0. Shifting of q_(gi) to s_(gi) has theadvantage that the process does not need to keep track of independentthresholds. Now the linear system of Eqn. (63) is solved as before. Ifthe solver leads to some (M) of the s_(gi,O) to be less than zero, theseare set equal to zero, and the remaining N−M−1 equations are solved. Theprocess is continued until optimal points are reached with alls_(gi,O)≥0.

FIGS. 5A-5E, collectively, is a flow chart that illustrates a controlscheme for gas lift operations for multiple wells based on directconstrained optimization using shifted gas flow rates. The operationscan be performed by a controller, such as the controller 162 of FIG. 3.

The operations begin in block 1001 where initial values for variablesand constraints N, Q_(g), θ_(gi) are provided by user input.

In block 1003, a shifted gas flow rate s_(gi) is defined for each well iwith respect to a gas injection flow rate q_(gi) and threshold flow rateθ_(gi) as: s_(gi)=q_(gi)−θ_(gi).

In block 1005, a constraint is defined based on the total allowable gasinjection flow rate Q_(g) and threshold gas injection flow rates θ_(gi)for each well i as:

${\overset{\sim}{Q}}_{g} = {Q_{g} - {\sum\limits_{i = 1}^{N}{\theta_{gi}.}}}$

In block 1007, an initial quadratic function {tilde over (q)}_(oi)(s_(gi)) is defined for each well i. The quadratic function {tilde over(q)}_(oi)(s_(gi)) specifies the oil production flow rate q_(oi) as afunction of s_(gi) for well i.

In block 1009, the control valve and injection line flow meter for welli (with i initially set to 1) is used to set two different q_(gi), whichcorrespond to two different s_(gi). For each one of the two q_(gi) (andcorresponding s_(gi)), the production line flow meter for the well i isused to measure the produced oil flow rate q_(oi) for non-monotonicq_(oi).

In block 1011, a predefined s_(gi), q_(oi) data point of (0,0) as wellas the two s_(gi), q_(oi) data points from block 1009 are used toconstruct the quadratic function {tilde over (q)}_(oi)(s_(gi)) for thewell i.

In block 1013, the quadratic function {tilde over (q)}_(oi)(s_(gi)) forthe well i of block 1011 is used to obtain an optimal shifted gas flowrate s_(qi,O) for well i, which assumes that no limits have been appliedby constraints. In this block, the optimal shifted gas flow rates_(qi,O) can be given by the coefficients of the quadratic function{tilde over (q)}_(oi)(s_(gi)) of Eqn. (57) as

$s_{{qi},O} = {\left( {- \frac{{\overset{\sim}{B}}_{i}}{2{\overset{\sim}{A}}_{i}}} \right).}$

In block 1015, the well index i is incremented.

In block 1017, the operations check that the current well index i doesnot exceed the maximum number i_(max) (which corresponds to the numberof N wells specified by user input in block 1001). If so, the operationsrevert to block 1007 to repeat the operations of blocks 1007 to 1015 forthe next well. If not, the operations continue to block 1019.

In block 1019, the operations check whether the condition

${\sum\limits_{i = 1}^{N}\left( s_{{gi},O} \right)} > {\overset{\sim}{Q}}_{g}$

is true or false, which checks that the sum of the optimal shifted gasflow rate s_(gi,O) for all wells does exceed {tilde over (Q)}_(g) (whichis based on the total available gas flow rate Q_(g) per Eqn. (61)). Iftrue, the operations continue to blocks 1021 to 1027 of FIGS. 5C and 5D.If false, the operations continue to blocks 1033 to 1041 of FIG. 5E.

In block 1021, the operations construct a computational model (based onEqn. (63)) for all N wells.

In block 1023, the operations solve the computational model of block1021 to determine the optimal shifted gas flow rate s_(gi,O) for the Nwells.

In block 1025, the operations identify a set of M wells where theoptimal shifted gas flow rate s_(gi,O) as determined in block 1023 isless than κ, and sets the shifted gas flow rate s_(gi) for these M wellsto 0.

In block 1027, the operations check whether M is greater than 0 (i.e.,the set of M wells is empty). If not (i.e., the set of M wells is notempty), the operations continue to block 1029 in order to reconstructthe computational model (Eqn. (63)) for N−M−1 wells with shifted gasflow rates s_(gi) specified according to block 1025 and the operationsrevert back to repeat the operations of blocks 1023 to 1027 for theupdated model. If so (i.e., the set of M wells is empty), the operationscontinue to block 1031.

In block 1031, the process has converged and the optimal shifted gasflow rate s_(gi,O) for all of the N wells is known, and the operationscontinue to blocks 1033 to 1041 of FIG. 5E.

In block 1033, the control valve and the injection line flow meter forwell i (with i initially set to 1) is used to set q_(gi) according tothe known s_(gi,O) (Eqn. (52)).

In block 1035, the production line flow meter for well i is used tomeasure q_(oi) for well i.

In block 1037, one of the three data points s_(gi), q_(oi) that was lastused to build the quadratic function {tilde over (q)}_(oi)(s_(gi)) forthe well i is dropped and replaced by the data point for s_(gi,O) andq_(oi) (as measured in block 1035) for the well i.

In block 1039, the well index i is incremented.

In block 1041, the operations check that the current well index i doesnot exceed the maximum number i_(max) (which corresponds to the numberof N wells specified by user input in block 1001). If so, the operationsrevert to block 1033 to repeat the operations of blocks 1033 to 1041 forthe next well. If not, the operations continue to block 1043 to wait fora predetermined wait time and then return to block 1007 to repeat thecontrol process. Note that the wait time of block 1043 is anticipated todepend on the surveillance rate suitable for continuous operation.

Nonlinear Optimization

In the foregoing, an efficient procedure to establish the relationshipof oil production with gas injection in a physical well wasdemonstrated. In addition, as the multiple wells in the gathering systemare treated in a separable manner, the outcome is that a definitivequadratic model can be realized for each concurrently over a number ofsampling steps. However, at this juncture, instead of using themulti-well scheme described above, one may resort to a broader and morerobust solution procedure enabled by application of nonlinearoptimization using a sequential quadratic programming (SQP) solver.Notably, under the stipulated assumptions, the models are smooth, convexand differentiable, and hence. the problem is simple to solve.

With regard to the above mentioned, consider a problem of N wells, eachwith a definitive realized production curve as per Eqn. (14) or Eqn.(57). For a given maximum available injection gas flow rate Q_(g), thegeneral nonlinear constrained optimization problem can be stated as

$\begin{matrix}{{{{\max\left( {\sum\limits_{i = 1}^{N}{q_{oi}\left( q_{gi} \right)}} \right)}\text{:}q_{gi}^{L}} \leq q_{gi} \leq q_{gi}^{U}},} & {{Eqn}.\mspace{14mu} (64)} \\{{{{subject}\mspace{14mu} {to}\mspace{14mu} {\sum\limits_{i = 1}^{N}q_{gi}}} \leq Q_{g}},} & {{Eqn}.\mspace{14mu} (65)}\end{matrix}$

where q_(gi) ^(L) and q_(gi) ^(U) indicate the lower and upper boundsfor gas injection in the i-th well, respectively. Notably, q_(gi) ^(L)may be zero, or set to θ_(gi) if a minimum injection is stipulated.q_(gi) ^(U) may be limited to the anticipated unconstrained optimum ofthe quadratic form as per Eqn. (8) since it will never be desirable toover inject. Indeed. if Σ_(i=1) ^(N)q_(gi,O) is less than the availablegas Q_(g), then the optimal solution is known directly from

$s_{{gi},O} = {- {\frac{{\overset{\sim}{B}}_{i}}{2{\overset{\sim}{A}}_{i}}.}}$

The advantage of solving the constrained nonlinear optimization problem(in place of the iterative procedure described above) is that itimplicitly provides the best possible solution subject to all stipulatedconstraints in one go (i.e., it will return the highest objectivemeasure given gas injection bounds and specified resource limitations).Moreover, any objective function can be considered, and any number andtype of additional constraints can be imposed. For example, well limitsby fluid phase, or collectively as cumulative quantities may bespecified. Hence, this procedure can impart greater flexibility,robustness and guarantees solution optimality for the problem posed.Finally, once established, the solution can be evaluated in practice,and the resulting sample points used to update the proxy quadraticmodels of each well in the field. The procedure can thus be repeated, asoutlined in the main section. The procedure however can employ acomputing environment (e.g., computer processing system) that isconfigured to process the SQP algorithm.

FIG. 6 illustrates adaptations to the control scheme of FIGS. 5A to 5Ebased on nonlinear optimization of the gas flow rates for multiplewells. The operations can be performed by a controller, such as thecontroller 162 of FIG. 3. In this control scheme, the iterativeoperations of blocks 1021 to 1029 of FIGS. 5C and 5D are replaced by theoperations 1021′ and 1023′ of FIG. 6.

In block 1021′, the operations construct a non-linear computation model(based on Eqns. (64) and (65)) for all N wells.

In block 1023′, the operations use the SQP solver to solve thenon-linear computation model of block 1023′ to determine the optimal gasflow rate q_(gi,O) for the N wells, and the operations continue to block1033 of FIG. 5E to update the gas injection rates for the for the Nwells based on the optimal gas flow rate q_(gi,O) for the N wells,update the quadratic function used for the optimization, and wait forpredetermined wait time before repeating the control process. Note thatthe use of shifted gas rates is unnecessary for the SQP solver and thusthe SQP solver can be configured to determine the optimal gas flow rateq_(gi,O) (non-shifted) for the N wells and the subsequent operations canuse the optimal gas flow rate q_(gi,O) for the N wells to update the gasinjection rates for the for the N wells.

Creep

A slow change in gas-lift performance is akin to creep, and is reflectedby the liquid production rate. By having a continuous update of thequadratic with changing oil production rate ensures a change in optimalgas injection rate, which in turn provides a data point at the ‘new’optimal. If the new rate allows us to update the quadratic by preservingnon-monotonicity, a continuous optimal operation of a creeping reservoiris possible.

Value Optimization

The gas lift control scheme that is described herein is configured tooptimize the flow rate of produced oil from one or more wells.Potentially more impactful, is the ability to carry out economicoptimization. Rather than optimizing q_(o) with respect to q_(g), onemay consider optimizing q_(o)−αq_(g)−βq_(w), or maximizingΣ(q_(oi)−αq_(gi)−βq_(wi)) (if Q_(g) limits gas), where a is a measure ofthe ratio of the well-head gas to liquid hydrocarbon value and has unitsof m³ ML⁻¹ and dimensionless β is the well-head water liability toliquid hydrocarbon value. Note that gas volume is at standardconditions, and that may show small variations from well to well becauseof access and location cost differences.

One advantage of value optimization is the direct inclusion of economicsin automation, without the need for elaborate simulation-based netpresent value calculation. The disadvantage is that the time value ofmoney is not considered, but this is balanced by not having to knowdetailed reservoir characteristics. The latter has large uncertainties.

The second advantage of value optimization is that a nearly flat q_(o)vs. q_(g) curve, close to optimal operating condition, ‘is sharpened’,enabling a more robust approach to the optimum. It also shifts theoptimal q_(g) to the left, and therefore is of value when Q_(g) limitstotal available gas.

Summary

The present disclosure describes a fully automated control system foroptimizing gas lift in one or more wells, which can include one or moreof the following operations:

-   -   Operate each well at a set of three gas injection flow rates        preferably ones for which the behavior of produced oil flow rate        vs. injected gas flow rate for the given well is non-monotonic.        Optionally, one of the three rates can be assumed as a produced        oil flow rate of zero at an injected gas flow rate of zero.    -   For a single well, the optimal gas injection flow rate is set        based on Eqn. (8), unless bounds are violated. In order to avoid        violating bounds it is preferable to employ data points of        produced oil flow rate vs. shifted gas flow rate (s_(g)) rather        than injected gas flow rate (q_(g)).    -   Use the new setting to update the quadratic function, preferably        retaining the non-monotonicity of q_(o) vs. q_(g) or s_(g).    -   If an anchor point well away from the set q_(g) (or s_(g)) is        detected, then adjust the set q_(g) (or s_(g)) towards the        anchor point by about 10% or less.    -   Upon converging to an optimal operation, periodically adjust the        flow rate, refining the quadratic function (which can involve        re-evaluating Eqn. (8)).    -   For multi-well control and optimization, when Q_(g) is not a        constraint, optimize each well based on s_(g) with

$s_{{gi},O} = {\left( {- \frac{{\overset{\sim}{B}}_{i}}{2{\overset{\sim}{A}}_{i}}} \right).}$

-   -   For multi-well control optimization, when constraints are        imposed, first solve unconstrained optimal settings, i.e.,

$s_{{gi},O} = {\left( {- \frac{{\overset{\sim}{B}}_{i}}{2{\overset{\sim}{A}}_{i}}} \right).}$

If

${{\sum\limits_{i = 1}^{N}\left( s_{{gi},O} \right)} \leq {\overset{\sim}{Q}}_{g}},$

then accept.

-   -   If

${{\sum\limits_{i = 1}^{N}\left( s_{{gi},O} \right)} > {\overset{\sim}{Q}}_{g}},$

impose {tilde over (Q)}_(g) constraint and re-evaluate the optimalsettings s_(gi,O) as per Eqn. (63).

-   -   If some (M) s_(gi,O)<0, then set them equal to zero, and        re-evaluate the optimal settings s_(gi,O) for N−M−1 equations.    -   Continue with new estimates of A_(i), B_(i) and C_(i) as in the        single well case; update the optimal settings until q_(gi) or        s_(gi) converge. Alternatively, use SQP with the updated        quadratics.    -   Reinitiate the entire process after periodically setting q_(gi)        to new values and reevaluating A_(i), B_(i) and C_(i), while        satisfying the Q_(g) constraint, and keeping new settings above        θ_(g) constraints.

FIG. 7 illustrates an example device 2500, with a processor 2502 andmemory 2504 that can be configured to implement various embodiments ofthe optimized gas lift control processes as discussed in thisdisclosure. For example, device 2500 can be configured to embody thecontroller 162 of FIGS. 2 and 3. Memory 2504 can also host one or moredatabases and can include one or more forms of volatile data storagemedia such as random-access memory (RAM), and/or one or more forms ofnonvolatile storage media (such as read-only memory (ROM), flash memory,and so forth).

Device 2500 is one example of a computing device or programmable deviceand is not intended to suggest any limitation as to scope of use orfunctionality of device 2500 and/or its possible architectures. Forexample, device 2500 can comprise one or more computing devices,programmable logic controllers (PLCs), etc.

Further, device 2500 should not be interpreted as having any dependencyrelating to one or a combination of components illustrated in device2500. For example, device 2500 may include one or more of computers,such as a laptop computer, a desktop computer, a mainframe computer,etc., or any combination or accumulation thereof.

Device 2500 can also include a bus 2508 configured to allow variouscomponents and devices, such as processors 2502, memory 2504, and localdata storage 2510, among other components, to communicate with eachother.

Bus 2508 can include one or more of any of several types of busstructures, including a memory bus or memory controller, a peripheralbus, an accelerated graphics port, and a processor or local bus usingany of a variety of bus architectures. Bus 2508 can also include wiredand/or wireless buses.

Local data storage 2510 can include fixed media (e.g., RAM, ROM, a fixedhard drive, etc.) as well as removable media (e.g., a flash memorydrive, a removable hard drive, optical disks, magnetic disks, and soforth).

One or more input/output (I/O) device(s) 2512 may also communicate via auser interface (UI) controller 2514, which may connect with I/Odevice(s) 2512 either directly or through bus 2508.

In one possible implementation, a network interface 2516 may communicateoutside of device 2500 via a connected network.

A media drive/interface 2518 can accept removable tangible media 2520,such as flash drives, optical disks, removable hard drives, softwareproducts, etc. In one possible implementation, logic, computinginstructions, and/or software programs comprising elements of module2506 may reside on removable media 2520 readable by mediadrive/interface 2518. Various processes of the present disclosure orparts thereof can be implemented by instructions and/or softwareprograms that are elements of module 2506. Such instructions and/orsoftware programs may reside on removable media 2520 readable by mediadrive/interface 2518 as is well known in the computing arts.

In one possible embodiment, input/output device(s) 2512 can allow a user(such as a human annotator) to enter commands and information to device2500, and also allow information to be presented to the user and/orother components or devices. Examples of input device(s) 2512 include,for example, sensors, a keyboard, a cursor control device (e.g., amouse), a microphone, a scanner, and any other input devices known inthe art. Examples of output devices include a display device (e.g., amonitor or projector), speakers, a printer, a network card, and so on.

Various processes of the present disclosure may be described herein inthe general context of software or program modules, or the techniquesand modules may be implemented in pure computing hardware. Softwaregenerally includes routines, programs, objects, components, datastructures, and so forth that perform particular tasks or implementparticular abstract data types. An implementation of these modules andtechniques may be stored on or transmitted across some form of tangiblecomputer-readable media. Computer-readable media can be any availabledata storage medium or media that is tangible and can be accessed by acomputing device. Computer readable media may thus comprise computerstorage media. “Computer storage media” designates tangible media, andincludes volatile and nonvolatile, removable and non-removable tangiblemedia implemented for storage of information such as computer readableinstructions, data structures, program modules, or other data. Computerstorage media include, but are not limited to, RAM, ROM, EEPROM, flashmemory or other memory technology, CD-ROM, digital versatile disks (DVD)or other optical storage, magnetic cassettes, magnetic tape, magneticdisk storage or other magnetic storage devices, or any other tangiblemedium which can be used to store the desired information, and which canbe accessed by a computer. Some of the methods and processes describedabove, can be performed by a processor. The term “processor” should notbe construed to limit the embodiments disclosed herein to any particulardevice type or system. The processor may include a computer system. Thecomputer system may also include a computer processor (e.g., amicroprocessor, microcontroller, digital signal processor, orgeneral-purpose computer) for executing any of the methods and processesdescribed above.

Some of the methods and processes described above, can be implemented ascomputer program logic for use with the computer processor. The computerprogram logic may be embodied in various forms, including a source codeform or a computer executable form. Source code may include a series ofcomputer program instructions in a variety of programming languages(e.g., an object code, an assembly language, or a high-level languagesuch as C, C++, or JAVA). Such computer instructions can be stored in anon-transitory computer readable medium (e.g., memory) and executed bythe computer processor. The computer instructions may be distributed inany form as a removable storage medium with accompanying printed orelectronic documentation (e.g., shrink wrapped software), preloaded witha computer system (e.g., on system ROM or fixed disk), or distributedfrom a server or electronic bulletin board over a communication system(e.g., the Internet or World Wide Web).

Alternatively or additionally, the processor may include discreteelectronic components coupled to a printed circuit board, integratedcircuitry (e.g., Application Specific Integrated Circuits (ASIC)),and/or programmable logic devices (e.g., a Field Programmable GateArrays (FPGA)). Any of the methods and processes described above can beimplemented using such logic devices.

There have been described and illustrated herein several embodiments ofprocesses and systems for controlling gas lift operations. Whileparticular embodiments have been described, it is not intended that theinvention be limited thereto, as it is intended that the invention be asbroad in scope as the art will allow and that the specification be readlikewise. It will therefore be appreciated by those skilled in the artthat yet other modifications could be made to the provided inventionwithout deviating from its spirit and scope as claimed.

What is claimed is:
 1. A method for controlling gas lift in at least oneproduction well, comprising: i) determining a quadratic function thatrelates oil flow rate produced from a production well to gas flow ratefor gas injected into the production well based on flow ratemeasurements for gas injected into the production well and correspondingflow rate measurements of oil produced from the production well; ii)using the quadratic function to determine an optimal flow rate for gasinjected into the production well; and iii) using the optimal flow ratefor gas injected into the production well to control the flow rate ofgas injected into the production well, and measuring an oil flow rateproduced from the production well at the optimal flow rate of gasinjected into the production well.
 2. A method according to claim 1,further comprising: repeating the operations of i) to iii), wherein forone or more iterations of i) to iii), the quadratic function of i) isbased on the measurement of oil flow rate in the previous iteration ofiii).
 3. A method according to claim 2, wherein: the measurement of oilflow rate in the previous iteration of iii) replaces a data point usedto determine the quadratic function in the previous iteration of iii)and preserves non-monotonicity of the quadratic function.
 4. A methodaccording to claim 2, wherein: the operations of i) to iii) are repeatedfor one or more iterations until a predetermined criterion is satisfied.5. A method according to claim 1, wherein: the quadratic functionrelates oil flow rate produced from the production well to injected gasflow rate for gas injected into the production well.
 6. A methodaccording to claim 1, wherein: the quadratic function relates oil flowrate produced from the production well to shifted gas flow rate for gasinjected into the production well, wherein shifted gas flow rate isbased on injected gas flow rate less a threshold gas flow rate.
 7. Amethod according to claim 1, wherein: the optimal flow rate for gasinjected into the production well as determined in ii) is based oncoefficients of the quadratic function.
 8. A method according to claim1, wherein: the operations of i) are performed separately for aplurality of production wells to determine a corresponding plurality ofquadratic functions, wherein each quadratic function relates oil flowrate produced from one of the plurality of production wells to gas flowrate for gas injected into that one production well; the operations ofii) use the plurality of quadratic functions to determine optimal flowrates for gas injected into each one of the plurality of productionwells; and the operations of iii) use the plurality of optimal flowrates for gas injected into plurality of the production wells to controlthe flow rate of gas injected into the plurality of production wells,and measure oil flow rates produced from the plurality of productionwells at the optimal flow rates of gas injected into the plurality ofproduction wells.
 9. A method according to claim 8, wherein: the optimalflow rate for gas injected into at least one of the plurality ofproduction wells as determined in ii) is based on coefficients of thequadratic function for the corresponding production well.
 10. A methodaccording to claim 8, wherein: the optimal flow rate for gas injectedinto at least one of the plurality of production wells as determined inii) is based on a constraint.
 11. A method according to claim 10,wherein: the constraint is based on a total available flow rate ofinjected gas from a common gas source shared by the plurality ofproduction wells.
 12. A method according to claim 10, wherein: theconstraint is based on a threshold rate of injected gas for at least oneproduction well.
 13. A method according to claim 10, wherein: theoptimal flow rates for gas injected into the plurality of productionwells as determined in ii) is based on a system of nonlinear equationssolved by a sequential quadratic programming (SQP) solver.
 14. A methodaccording to claim 1, wherein: the quadratic function of i) relates oilflow rate produced from the production well to gas flow rate for gasinjected into the production well at standard conditions.
 15. A systemfor controlling gas lift in at least one production well, comprising: afirst flow meter for measuring flow rate of gas injected into aproduction well; a control valve for controlling flow rate of gasinjected into the production well; a second flow meter for measuringflow rate of oil produced from the production well; and a controller,operably coupled to the first flow meter, the control valve, and thesecond flow meter, wherein the controller is configured to: i) determinea quadratic function that relates oil flow rate produced from aproduction well to gas flow rate for gas injected into the productionwell based on flow rate measurements for gas injected into theproduction well performed by the first flow meter and corresponding flowrate measurements of oil produced from the production well performed bythe second flow meter; ii) use the quadratic function to determine anoptimal flow rate for gas injected into the production well; and iii)use the optimal flow rate for gas injected into the production well tocontrol the control valve to control the flow rate of gas injected intothe production well, and obtain a measurement of oil flow rate producedfrom the production well at the optimal flow rate of gas injected intothe production well performed by the second flow meter.
 16. A systemaccording to claim 15, wherein: the controller is further configured torepeat the operations of i) to iii), wherein for one or more iterationsof i) to iii), the quadratic function of i) is based on the measurementof oil flow rate performed by the second flow meter in the previousiteration of iii).
 17. A system according to claim 16, wherein: thecontroller is further configured such that the measurement of oil flowrate in the previous iteration of iii) replaces a data point used todetermine the quadratic function in the previous iteration of iii) andpreserves non-monotonicity of the quadratic function.
 18. A systemaccording to claim 16, wherein: the controller is further configured torepeat the operations of i) to iii) for one or more iterations until apredetermined criterion is satisfied.
 19. A system according to claim15, wherein: the quadratic function relates oil flow rate produced fromthe production well to injected gas flow rate for gas injected into theproduction well.
 20. A system according to claim 15, wherein: thequadratic function relates oil flow rate produced from the productionwell to shifted gas flow rate for gas injected into the production well,wherein shifted gas flow rate is based on injected gas flow rate less athreshold gas flow rate.
 21. A system according to claim 15, wherein:the controller is further configured such that the optimal flow rate forgas injected into the production well as determined in ii) is based oncoefficients of the quadratic function.
 22. A system according to claim15, wherein: the first flow meter, the control valve, and the secondflow meter are provided separately for a plurality of production wells;the controller is operably coupled to the first flow meter, the controlvalve, and the second flow meter for the plurality of production wells;and the controller is further configured such that: the operations of i)are performed separately for the plurality of production wells todetermine a corresponding plurality of quadratic functions, wherein eachquadratic function relates oil flow rate produced from one of theplurality of production wells to gas flow rate for gas injected intothat one production well, the operations of ii) use the plurality ofquadratic functions to determine optimal flow rates for gas injectedinto each one of the plurality of production wells, and the operationsof iii) use the plurality of optimal flow rates for gas injected intoplurality of the production wells to control the control valve tocontrol flow rate of gas injected into the plurality of productionwells, and measure oil flow rates produced from the plurality ofproduction wells at the optimal flow rates of gas injected into theplurality of production wells.
 23. A system according to claim 22,wherein: the optimal flow rate for gas injected into at least one of theplurality of production wells as determined in ii) is based oncoefficients of the quadratic function for the corresponding productionwell.
 24. A system according to claim 22, wherein: the optimal flow ratefor gas injected into at least one of the plurality of production wellsas determined in ii) is based on a constraint.
 25. A system according toclaim 24, wherein: the constraint is based on a total available flowrate of injected gas from a common gas source shared by the plurality ofproduction wells.
 26. A system according to claim 24, wherein: theconstraint is based on a threshold rate of injected gas for at least oneproduction well.
 27. A system according to claim 24, wherein: theoptimal flow rates for gas injected into the plurality of productionwells as determined in ii) is based on a system of nonlinear equationssolved by a sequential quadratic programming (SQP) solver.
 28. A systemaccording to claim 15, wherein: the controller comprises a processor.29. A method for controlling gas lift in at least one production well,comprising: i) determining a quadratic function that relates value ofoil produced from a production well to gas flow rate for gas injectedinto the production well based on flow rate measurements for gasinjected into the production well and corresponding flow ratemeasurements of oil produced from the production well; ii) using thequadratic function to determine an optimal flow rate for gas injectedinto the production well; and iii) using the optimal flow rate for gasinjected into the production well to control the flow rate of gasinjected into the production well, and measuring an oil flow rateproduced from the production well at the optimal flow rate of gasinjected into the production well.
 30. A method according to claim 29,further comprising: repeating the operations of i) to iii), wherein forone or more iterations of i) to iii), the quadratic function of i) isbased on the measurement of oil flow rate in the previous iteration ofiii).
 31. A method according to claim 30, wherein: the measurement ofoil flow rate in the previous iteration of iii) is used to replace adata point used to determine the quadratic function in the previousiteration of iii) and preserves non-monotonicity of the quadraticfunction.
 32. A method according to claim 30, wherein: the operations ofi) to iii) are repeated for one or more iterations until a predeterminedcriterion is satisfied.